import numpy as np
import pdb
from scipy.spatial import Delaunay, ConvexHull
class room_mesh_generate():
    def  __init__(self,new_points):
        super().__init__()
        self.points = new_points
        self.boundary = self.boundary_identification()
        # in and out 
        self.inlet = np.array(self.inlet_identification())
        self.outlet = np.array(self.outlet_identification())
        self.inlet_bd = np.array(self.inlet_bd_identification())
        self.outlet_bd = np.array(self.outlet_bd_identification())
        # pdb.set_trace()
        idx_in_bd = np.where(np.all(np.isin(self.inlet,self.inlet_bd),axis=1))[0]
        idx_out_bd = np.where(np.all(np.isin(self.outlet,self.outlet_bd),axis=1))[0]
        self.inlet_inn = np.delete(self.inlet,idx_in_bd,axis=0)
        self.outlet_inn = np.delete(self.outlet,idx_out_bd,axis=0)
        idx_inlet_inn = np.where(np.all(np.isin(self.boundary,self.inlet_inn),axis=1))[0]
        idx_outlet_inn = np.where(np.all(np.isin(self.boundary,self.outlet_inn),axis=1))[0]
        self.wall = np.delete(self.boundary,np.concatenate([idx_inlet_inn, idx_outlet_inn]),axis=0)
        self.internal = np.delete(self.points,np.where(np.all(np.isin(self.points,self.boundary),axis=1))[0],axis=0)
        self.points = np.concatenate([self.internal,self.inlet_inn,self.outlet_inn,self.wall],axis=0)
        self.idx_Inlet,self.idx_Outlet,self.idx_Wall =self.extract_id()
        self.Tri = Delaunay(self.points,qhull_options='QJ')
        self.faces = self.extract_triangle()
        print('节点数量:',self.points.shape[0])
        print('四面体数量:',self.Tri.simplices.shape[0])
        print('面单元数量:',self.faces.shape[0])
    
    # extract faces
    def extract_triangle(self):
        simplices = self.Tri.simplices
        triangles = []
        for i in range(simplices.shape[0]):
            idx = simplices[i]
            trianges1 = [idx[0], idx[1], idx[2]]
            trianges2 = [idx[0], idx[1], idx[3]]
            trianges3 = [idx[0], idx[2], idx[3]]
            trianges4 = [idx[1], idx[2], idx[3]]
            trianges1.sort()
            trianges2.sort()
            trianges3.sort()
            trianges4.sort()
            triangles.append([trianges1,trianges2,trianges3,trianges4])
        triangles = np.asarray(triangles).reshape(-1,3)
        triangles = np.unique(triangles,axis = 0)
        return  triangles
    # extract boundary
    def boundary_identification(self):
        points = self.points
        boundary = []
        for i in range(points.shape[0]):
            p=points[i].tolist()
            a = (1.4999<=p[0]<=1.50001) or (-1.50001<=p[0]<=-1.499999999)
            b = (1.4999<=p[1]<=1.50001) or (-1.50001<=p[1]<=-1.499999999)
            c = (1.4999<=p[2]<=1.50001) or (-1.50001<=p[2]<=-1.499999999)
            if (a or b or c):
                boundary.append(p)        
        return np.array(boundary)
    #  extract outlet
    def outlet_identification(self):
        points = self.points
        outlet = []
        for i in range(points.shape[0]):
            p = points[i].tolist()
            if (1.4999999<=p[1]<=1.5000001):
                r = np.linalg.norm(points[i]-np.array([0.0,1.5,0.0]))
                if r <= 0.100001:
                    outlet.append(p)
        return outlet
    #  extract inlet
    def inlet_identification(self):
        points = self.points
        outlet = []
        for i in range(points.shape[0]):
            p = points[i].tolist()
            if (-1.5000001<=p[0]<=-1.4999999):
                if (0.9499999<=p[1]<=1.05000001):
                    if (-0.499999<=p[2]<=0.5000001):
                        outlet.append(p)
        return outlet
    #  extract outlet boundary
    def outlet_bd_identification(self):
        points = self.points
        outlet = []
        for i in range(points.shape[0]):
            p = points[i].tolist()
            if (1.4999999<=p[1]<=1.5000001):
                r = np.linalg.norm(points[i]-np.array([0.0,1.5,0.0]))
                if 0.0999999< r < 0.100001:
                    outlet.append(p)
        return outlet
    #  extract inlet boundary
    def inlet_bd_identification(self):
        points = self.points
        inlet = []
        for i in range(points.shape[0]):
            p = points[i].tolist()
            if (-1.5000001<p[0]<-1.4999999):
                a = ((0.9499999<p[1]<0.95000001) or (1.0499999<p[1]<1.05000001))
                b = ((-0.5000001<p[2]<-0.499999) or (0.499999<p[2]<0.5000001))
                if a or b:
                    inlet.append(p)
        return inlet
    
    #  extract points_set index
    def extracting(self,data_all,indexing_data):
        indxing = []
        for i in range(indexing_data.shape[0]):
            s = data_all-indexing_data[i]
            judge = np.count_nonzero(s,axis=1)
            idx = (np.where(judge==0)[0])
            idx = idx.tolist()
            indxing +=idx
        return indxing
    
    #  extract index
    def extract_id(self):
        # wall
        idx_Inlet = self.extracting(self.points,self.inlet)
        idx_Outlet = self.extracting(self.points,self.outlet)
        idx_Wall = self.extracting(self.points,self.wall)
        return idx_Inlet,idx_Outlet,idx_Wall
    
    #  append faces information
    def append_left_right(self):
        faces = self.faces
        simplices = self.Tri.simplices
        info_faces = np.empty(shape=(faces.shape[0],5),dtype = 'int64')
        for i in range(faces.shape[0]):
            f = faces[i]
            judge = np.count_nonzero(np.isin(simplices,faces[i]),axis=1)
            s = np.where(judge == 3)[0]
            for j in s:
                idx = np.where(np.isin(simplices[j],faces[i],invert=True))[0]
                left = simplices[j][idx]
                f = np.append(f,left)
            if f.shape[0]==4:
                f = np.append(f,np.ones(1)*-1)
            info_faces[i] = f
            if i%10000 ==0:
                print('已完成faces计算的数量：',i/faces.shape[0]*100,'%')
        self.info_faces = info_faces + 1
        return self.info_faces 

    #  sorting the faces information
    def sorting_faces(self):
        faces = self.info_faces
        idx_Inlet,idx_Outlet,idx_Wall= self.extract_id()
        Inlet_faces = []
        Outlet_faces = []
        Wall_faces = []
        Internal_faces = []
        for i in range(faces.shape[0]):
            Inlet_judge = np.all(np.isin((faces[i,0:3]-1),idx_Inlet))
            Outlet_judge = np.all(np.isin((faces[i,0:3]-1),idx_Outlet))
            Wall_judge = np.all(np.isin((faces[i,0:3]-1),idx_Wall))
            if Inlet_judge:
                Inlet_faces.append(faces[i])
            elif Outlet_judge:
                Outlet_faces.append(faces[i])
            elif Wall_judge:
                Wall_faces.append(faces[i])
            else:
                Internal_faces.append(faces[i])
            if i%10000 ==0:
                print('已完成faces计算的数量：',i/faces.shape[0]*100,'%')
        self.Inlet_faces = Inlet_faces
        self.Outlet_faces = Outlet_faces
        self.Wall_faces = Wall_faces
        self.Internal_faces = Internal_faces
        
        return Inlet_faces,Outlet_faces,Wall_faces,Internal_faces
    
    def write_msh(self,dim = 3,facetype =3):
        Tri_All = self.Tri
        points = self.points
        idx_Inlet = self.idx_Inlet
        idx_Outlet = self.idx_Outlet
        idx_Wall = self.idx_Wall
        Internal_faces = np.array(self.Internal_faces,dtype='int64')
        Inlet_faces = np.array(self.Inlet_faces,dtype='int64')
        Outlet_faces = np.array(self.Outlet_faces,dtype='int64')
        Wall_faces = np.array(self.Wall_faces,dtype='int64')
        info_faces = np.concatenate([Internal_faces,Inlet_faces,Outlet_faces,Wall_faces],axis = 0)

        title = ['(0 grid written by ANSYS Meshing \n',
                 'nodes:       (10 (id start end type) (x y z ...))\n',
                 'faces:       (13 (id start end type elemType)\n',
                 '   (v-0 v-1 .. v-n right-cell left-cell ...))\n',
                 'cells:       (12 (id start end type elemtype))\n',
                 'parent-face: (59 (start end parent child) (nchilds child0 child1 ...))\n',
                 ')\n']

        num_simplices = Tri_All.simplices.shape[0]
        num_Points = points.shape[0]
        num_faces = info_faces.shape[0]
        num_face_inlet = np.array(Inlet_faces).shape[0]
        num_face_outlet = np.array(Outlet_faces).shape[0]
        num_face_wall = np.array(Wall_faces).shape[0]
        num_face_internal = np.array(Internal_faces).shape[0]


        id_Nodes_bd = np.concatenate([idx_Inlet,idx_Outlet,idx_Wall])
        points_boundary = points[id_Nodes_bd,:]
        points_internal = np.delete(points,id_Nodes_bd,axis=0)
        num_internal = points_internal.shape[0]
        num_booundary = points_boundary.shape[0]

        zoneid_inn = 1
        zoneid_inlet = 5
        zoneid_outlet = 6
        zoneid_wall = 7
        
        with open('room_mesh.msh', 'w') as f:
            # write title
            for i in title:
                f.write(i)
            # write dim
            f.write('(2 ' + str(dim) +')\n')
            # write basic info
            ## Nodes info
            a = hex(num_Points)[2:]
            f.write('(10 (0 1 '+ str(a) +' 0))\n')
            ## Faces info
            a = hex(num_faces)[2:]
            f.write('(13 (0 1 '+ str(a) +' 0))\n')
            ## Cells info
            a = hex(num_simplices)[2:]
            f.write('(12 (0 1 '+ str(a) +' 0))\n')

            # write Nodes message
            ## write internal nodes
            a = hex(num_internal)[2:]
            f.write('(10 (3 1 ' + str(a) + ' 1 ' +str(dim)+')(\n')
            for i in points_internal:
                f.write(str(i).replace('[', '').replace(']', '') +'\n')
            f.write('))\n')
            ## write boundary nodes
            a = hex(num_internal+1)[2:]
            b = hex(num_Points)[2:]
            f.write('(10 (3 ' + str(a) +' ' + str(b) + ' 1 ' + str(dim) + ')(\n')
            for i in points_boundary:
                f.write(str(i).replace('[', '').replace(']', '') +'\n')
            f.write('))\n')

            # write Faces message

            ## write internal faces
            a = hex(num_face_internal)[2:]
            f.write('(13 ( '+str(zoneid_inn)+' 1 ' + str(a) + ' 2 ' +str(facetype)+')(\n')
            for i in Internal_faces:
                hexval = []
                for j in i:
                    hexval.append(hex(j)[2:])
                # pdb.set_trace()
                f.write(str(hexval).replace('[', '').replace(']', '').replace('\'', '').replace(',', '') +'\n')
            f.write('))\n')

            ## write inlet faces
            a = hex(num_face_internal + 1)[2:]
            b = hex(num_face_internal + num_face_inlet + num_face_outlet)[2:]
            f.write('(13 ( '+str(zoneid_inlet)+ str(a) +' ' + str(b) + ' a ' +str(facetype)+')(\n')
            for i in Inlet_faces:
                hexval = []
                for j in i:
                    hexval.append(hex(j)[2:])
                # pdb.set_trace()
                f.write(str(hexval).replace('[', '').replace(']', '').replace('\'', '').replace(',', '') +'\n')
            f.write('))\n')

            ## write outlet faces
            a = hex(num_face_internal+num_face_inlet+1)[2:]
            b = hex(num_face_internal+num_face_inlet)[2:]
            f.write('(13 ( '+str(zoneid_outlet)+ str(a) +' ' + str(b) + ' 5 ' +str(facetype)+')(\n')
            for i in Outlet_faces:
                hexval = []
                for j in i:
                    hexval.append(hex(j)[2:])
                # pdb.set_trace()
                f.write(str(hexval).replace('[', '').replace(']', '').replace('\'', '').replace(',', '') +'\n')
            f.write('))\n')

            ## write Wall faces
            a = hex(num_face_internal+num_face_inlet+num_face_outlet+1)[2:]
            b = hex(num_face_internal+num_face_inlet+num_face_outlet+num_face_wall)[2:]
            f.write('(13 ( '+str(zoneid_wall)+ str(a) +' ' + str(b) + ' 5 ' +str(facetype)+')(\n')
            for i in Wall_faces:
                hexval = []
                for j in i:
                    hexval.append(hex(j)[2:])
                # pdb.set_trace()
                f.write(str(hexval).replace('[', '').replace(']', '').replace('\'', '').replace(',', '') +'\n')
            f.write('))\n')

            # write Cells 
            a = hex(num_simplices)[2:]
            f.write('(12 (2 1 '+ str(a) +' 1 2))\n')

            ## write zone info
            f.write('(45 (1 interior interior-fff___)())\n')
            f.write('(45 (2 fluid fff___)())\n')
            f.write('(45 (5 velocity-inlet inlet)())\n')
            f.write('(45 (6 pressure-outlet outlet)())\n')
            f.write('(45 (7 wall roof)())\n')